import numpy as np
def P_to_B(P_matrix, centering=False):
    P_x = np.sum(P_matrix, axis=1)
    P_y = np.sum(P_matrix, axis=0)
    if(centering==False):
        return np.diag(1/np.sqrt(P_x)) @ P_matrix @ np.diag(1/np.sqrt(P_y))
    dim = P_x.shape[0]
    offset_m = np.kron(P_x, P_y).reshape((dim, dim))
    return np.diag(1/np.sqrt(P_x)) @ (P_matrix - offset_m) @ np.diag(1/np.sqrt(P_y))
        
def mutual_info(P_matrix):
    P_x = np.sum(P_matrix, axis=1)
    P_y = np.sum(P_matrix, axis=0)
    dim = P_x.shape[0]
    offset_m = np.kron(P_x, P_y).reshape((dim, dim))
    return np.sum(P_matrix * np.log(P_matrix / offset_m))